{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Coords 2: Transforming between coordinate systems\n",
    "\n",
    "## Authors\n",
    "Erik Tollerud, Kelle Cruz, Stephen Pardy, Stephanie T. Douglas\n",
    "\n",
    "## Learning Goals\n",
    "- Create `astropy.coordinates.SkyCoord` objects\n",
    "- Transform to different coordinate systems on the sky\n",
    "- Transform to altitude/azimuth coordinates from a specific observing site\n",
    "\n",
    "## Keywords\n",
    "coordinates, units, observational astronomy\n",
    "\n",
    "## Summary\n",
    "\n",
    "In this tutorial we demonstrate how to define astronomical coordinates using the `astropy.coordinates` \"frame\" classes. We then show how to transform between the different built-in coordinate frames, such as from ICRS (RA, Dec) to Galactic (l, b). Finally, we show how to compute altitude and azimuth from a specific observing site."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Imports"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Third-party dependencies\n",
    "from astropy import units as u\n",
    "from astropy.coordinates import SkyCoord\n",
    "import numpy as np"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Set up matplotlib and use a nicer set of plot parameters\n",
    "from astropy.visualization import astropy_mpl_style\n",
    "import matplotlib.pyplot as plt\n",
    "plt.style.use(astropy_mpl_style)\n",
    "%matplotlib inline"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Section 0: Quickstart"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "<div class=\"alert alert-info\">\n",
    "\n",
    "**Note:** If you already worked through [Coords 1](http://learn.astropy.org/rst-tutorials/Coordinates-Intro.html?highlight=coordinates) you can feel free to skip to [Section 1](#Section-1:).\n",
    "\n",
    "</div>"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "In Astropy, the most common object you'll work with for coordinates is `SkyCoord`.  A `SkyCoord` can most easily be created directly from angles as shown below.  \n",
    "\n",
    "In this tutorial we'll be converting between frames. Let's start in the ICRS frame (which happens to be the default.) \n",
    "\n",
    "For much of this tutorial we'll work with the Hickson Compact Group 7. We can create an object either by passing the degrees explicitly (using the astropy [units](http://docs.astropy.org/en/stable/units/index.html) library) or by passing in strings. The two coordinates below are equivalent:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "hcg7_center = SkyCoord(9.81625*u.deg, 0.88806*u.deg, frame='icrs')  # using degrees directly\n",
    "print(hcg7_center)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "hcg7_center = SkyCoord('0h39m15.9s', '0d53m17.016s', frame='icrs')  # passing in string format\n",
    "print(hcg7_center)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "We can get the right ascension and declination components of the object directly by accessing those attributes."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "print(hcg7_center.ra)\n",
    "print(hcg7_center.dec)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Section 1: "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Introducing frame transformations"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "`astropy.coordinates` provides many tools to transform between different coordinate systems. For instance, we can use it to transform from ICRS coordinates (in RA and Dec) to Galactic coordinates.\n",
    "\n",
    "To understand the code in this section, it may help to read over the [overview of the astropy coordinates scheme](http://astropy.readthedocs.org/en/latest/coordinates/index.html#overview-of-astropy-coordinates-concepts).  The key piece to understand is that all coordinates in Astropy are in particular \"frames\" and we can transform between a specific `SkyCoord` object in one frame to another.  For example, we can transform our previously-defined center of HCG 7 from ICRS to Galactic coordinates:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "hcg7_center = SkyCoord(9.81625*u.deg, 0.88806*u.deg, frame='icrs')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "There are three different ways of transforming coordinates. Each has its pros and cons, but all should give you the same result. The first way to transform to other built-in frames is by specifying those attributes. For instance, let's see the location of HCG 7 in Galactic coordinates."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### Transforming coordinates using attributes:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "hcg7_center.galactic"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### Transforming coordinates using the transform_to() method and other coordinate object"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The above is actually a special \"quick-access\" form that internally does the same as what's in the cell below: it uses the [`transform_to()`](http://docs.astropy.org/en/stable/api/astropy.coordinates.SkyCoord.html#astropy.coordinates.SkyCoord.transform_to) method to convert from one frame to another. We can pass in an empty coordinate class to specify what coordinate system to transform into."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "from astropy.coordinates import Galactic  # new coordinate baseclass\n",
    "hcg7_center.transform_to(Galactic())"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### Transforming coordinates using the transform_to() method and a string"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Finally, we can transform using the `transform_to()` method and a string with the name of a built-in coordinate system."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "hcg7_center.transform_to('galactic')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "We can transform to many coordinate frames and equinoxes. \n",
    "\n",
    "These coordinates are available by default:\n",
    "\n",
    " - ICRS\n",
    " - FK5\n",
    " - FK4\n",
    " - FK4NoETerms\n",
    " - Galactic\n",
    " - Galactocentric\n",
    " - Supergalactic\n",
    " - AltAz\n",
    " - GCRS\n",
    " - CIRS\n",
    " - ITRS\n",
    " - HCRS\n",
    " - PrecessedGeocentric\n",
    " - GeocentricTrueEcliptic\n",
    " - BarycentricTrueEcliptic\n",
    " - HeliocentricTrueEcliptic\n",
    " - SkyOffsetFrame\n",
    " - GalacticLSR\n",
    " - LSR\n",
    " - BaseEclipticFrame\n",
    " - BaseRADecFrame\n",
    "\n",
    "Let's focus on just a few of these. We can try FK5 coordinates next:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "hcg7_center_fk5 = hcg7_center.transform_to('fk5')\n",
    "print(hcg7_center_fk5)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "And, as with the Galactic coordinates, we can acheive the same result by importing the FK5 class from the `astropy.coordinates` package. This also allows us to change the equinox."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "from astropy.coordinates import FK5\n",
    "hcg7_center_fk5.transform_to(FK5(equinox='J1975'))  # precess to a different equinox  "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "<div class=\"alert alert-warning\">\n",
    "\n",
    "**Beware:** Changing frames also changes some of the attributes of the object, but usually in a way that makes sense. The following code should fail.\n",
    "\n",
    "</div>"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "tags": [
     "raises-exception"
    ]
   },
   "outputs": [],
   "source": [
    "hcg7_center.galactic.ra  # should fail because Galactic coordinates are l/b not RA/Dec"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Instead, we now have access to the l and b attributes:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "print(hcg7_center.galactic.l, hcg7_center.galactic.b)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Section 2:"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Transform frames to get to altitude-azimuth (\"AltAz\")\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "To actually do anything with observability we need to convert to a frame local to an on-earth observer.  By far the most common choice is horizontal altitude-azimuth coordinates, or \"AltAz\".  We first need to specify both where and when we want to try to observe.\n",
    "\n",
    "We'll need to import a few more specific modules:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "from astropy.coordinates import EarthLocation\n",
    "from astropy.time import Time"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Let's first see the sky position at Kitt Peak National Observatory in Arizona."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Kitt Peak, Arizona\n",
    "kitt_peak = EarthLocation(lat='31d57.5m', lon='-111d35.8m', height=2096*u.m)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "For known observing sites we can enter the name directly."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "kitt_peak = EarthLocation.of_site('Kitt Peak')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "We can see the list of observing sites:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "EarthLocation.get_site_names()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Let's check the altitude at 1 AM UTC, which is 6 PM AZ mountain time:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "observing_time = Time('2010-12-21 1:00')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Now we use these to create an `AltAz` frame object.  Note that this frame has some other information about the atmosphere, which can be used to correct for atmospheric refraction.  Here we leave that alone, because the default is to ignore this effect (by setting the pressure to 0)."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "from astropy.coordinates import AltAz\n",
    "\n",
    "aa = AltAz(location=kitt_peak, obstime=observing_time)\n",
    "print(aa)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Now we can transform our ICRS `SkyCoord` to `AltAz` to get the location in the sky over Kitt Peak at the requested time."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "hcg7_center.transform_to(aa)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "To look at just the altitude we can `alt` attribute:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "hcg7_center.transform_to(aa).alt"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Alright, it's at 55 degrees at 6 PM, but that's pretty early to be observing.  We could try various times one at a time to see if the airmass is at a darker time, but we can do better: let's try to create an airmass plot."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# this gives a Time object with an *array* of times\n",
    "delta_hours = np.linspace(0, 6, 100)*u.hour\n",
    "full_night_times = observing_time + delta_hours\n",
    "full_night_aa_frames = AltAz(location=kitt_peak, obstime=full_night_times)\n",
    "full_night_aa_coos = hcg7_center.transform_to(full_night_aa_frames)\n",
    "\n",
    "plt.plot(delta_hours, full_night_aa_coos.secz)\n",
    "plt.xlabel('Hours from 6pm AZ time')\n",
    "plt.ylabel('Airmass [Sec(z)]')\n",
    "plt.ylim(0.9,3)\n",
    "plt.tight_layout()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Great!  Looks like the lowest airmass is in another hour or so (7 PM).  But that might still be twilight... When should we start observing for proper dark skies?  Fortunately, Astropy provides a ``get_sun`` function that can be used to check this.  Let's use it to check if we're in 18-degree twilight or not."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "from astropy.coordinates import get_sun\n",
    "\n",
    "full_night_sun_coos = get_sun(full_night_times).transform_to(full_night_aa_frames)\n",
    "plt.plot(delta_hours, full_night_sun_coos.alt.deg)\n",
    "plt.axhline(-18, color='k')\n",
    "plt.xlabel('Hours from 6pm AZ time')\n",
    "plt.ylabel('Sun altitude')\n",
    "plt.tight_layout()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Looks like it's just below 18 degrees at 7 PM, so you should be good to go!"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "We can also look at the object altitude at the present time and date:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "now = Time.now()\n",
    "hcg7_center = SkyCoord(9.81625*u.deg, 0.88806*u.deg, frame='icrs')\n",
    "kitt_peak_aa = AltAz(location=kitt_peak, obstime=now)\n",
    "print(hcg7_center.transform_to(kitt_peak_aa))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Exercises"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Exercise 1"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Try to compute to some arbitrary precision (rather than eyeballing on a plot) when 18 degree twilight or sunrise/sunset hits on that night."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Exercise 2"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Try converting the HCG 7 coordinates to an equatorial frame at some other equinox a while in the past (like J2000).  Do you see the precession of the equinoxes?\n",
    "\n",
    "Hint: To see a diagram of the supported frames look [here](http://docs.astropy.org/en/stable/coordinates/#module-astropy.coordinates) or on the list above. One of those will do what you need if you give it the right frame attributes."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Exercise 3"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Try looking at the altitude of HCG 7 at another observatory."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Wrap-up"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "For more documentation on the many other features of `astropy.coordinates`, check out [its section of the documentation](http://astropy.readthedocs.org/en/latest/coordinates/index.html).\n",
    "\n",
    "You might also be interested in [the astroplan affiliated package](http://astroplan.readthedocs.org/), which uses the `astropy.coordinates` to do more advanced versions of the tasks in the last section of this tutorial."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": []
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.6.7"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}
